Effect of disorder on the ground-state properties of graphene 
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We calculate the ground-state energy of Dirac electrons in graphene in the presence of disorder. We 
take randomly distributed charged impurities at a fixed distance from the graphene sheet and surface 
fluctuations (ripples) as the main scattering mechanisms. Mode-coupling approach to scattering rate 
and random-phase approximation for ground-state energy incorporating the many-body interactions 
and the disorder effects yields good agreement with experimental inverse compressibility. 
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I. INTRODUCTION 

Two-dimensional crystals of Carbon atoms (graphene) 
are recently discovered 1 . Graphene, a single, one-atom 
thick sheet of carbon atoms arranged in a honeycomb 
lattice. High quality graphene single crystals some thou- 
sands of fim 2 in size are sufficient for most fundamen- 
tal physics studies. 2 There are significant efforts to grow 
graphene epitaxially 3 - by thermal decomposition of SiC, 
or by vapor deposition of hydrocarbons on catalytic 
metallic surfaces which could later be etched away leav- 
ing graphene on an insulating substrate. 

This stable crystal has attracted considerable atten- 
tion because of its unusual effective many-body prop- 
erties 4 , quasi-particle properties and its Landau Fermi 
liquid pictured and the effect of electron-electron inter- 
actions to plasmon behavior and angle resolved photoe- 
mission spectroscopy (ARPES) 6 that follow from chiral 
band states and because of potential applications. The 
low energy quasi-particle excitations in graphene are lin- 
early dispersing, described by Dirac cones at the edges of 
the first Brillouin zone. It is very hard for alien atoms to 
replace the carbon atoms in the graphene structure be- 
cause of the robustness and specificity of the a bonding. 
Due to that, electron mean-free path I, in graphene can 
be very large. One of the important issues in graphene 
is its quantum transport properties having the universal 
minimum conductivity at the Dirac point. Initially, it 
was believed that this universality is a native property 7 
but recent experimental^ and theoretica l 10 ' 11 ' 12 ' 13 ' 14 re- 
ports indicate that the transport properties are very sen- 
sitive to impurities and defects and minimum conductiv- 
ity is not universal. 

Conventional two-dimensional electron gas (2DEG) 
has been a fertile source of surprising new physics for 
more than four decades. Although the exploration of 
graphene is still at an early stage, it is already clear— that 
the strong field properties of Dirac electrons in graphene 
are different from and as rich as those of a semiconduc- 
tor heterojunction 2DEG. The Fermi liquid phenomenol- 
ogy of Dirac electrons in graphene^ and conventional 
2DEG — have the same structure, since both systems are 



isotropic and have a single circular Fermi surface. The 
strength of interaction effects in a conventional 2DEG 
increases with decreasing carrier density. At low den- 
sities, the quasiparticle weight Z is small, the velocity 
is suppressed 1 ^, the charge compressibility changes sign 
from positive to negative 1 ^, and the spin-susceptibility is 
strongly enhanced 1 ^. These effects emerge from an inter- 
play between exchange interactions and quantum fluctu- 
ations of charge and spin in the 2DEG. 

In the Dirac electrons in graphene, it was shown 4 *^ 
that interaction effects also become noticeable with de- 
creasing density, although more slowly, that the quasi- 
particle weight Z tends to larger values, that the velocity 
is enhanced rather than suppressed, and that the influ- 
ence of interactions on the compressibility and the spin- 
susceptibility changes sign. These qualitative differences 
are due to exchange interactions between electrons near 
the Fermi surface and electrons in the negative energy sea 
and to interband contributions to Dirac electrons from 
charge and spin fluctuations. 

Compressibility measurements of conventional 2DEG 
have been carried ov^M- and it is found qualitatively that 
Coulomb interactions affect the compressibility at suffi- 
ciently low electron density or strong coupling constant 
region. Recently, the local compressibility of graphene 
has been measured 18 using a scannable single electron 
transistor and it is argued that the measured compress- 
ibility is well described by the kinetic energy contribution 
and it is suggested that exchange and correlation effects 
have canceling contributions. From the theoretical point 
of view, the compressibility was first calculated by Peres 
et al*^ considering the exchange contribution to the non- 
interacting doped or undoped graphene flake. A related 
quantity d/i/dn (where \i is the chemical potential and n 
is the electron density) is recently considered by Hwang 
et al£^ within the same approximation. Going beyond 
the exchange contribution, the correlation effects were 
taken into account by Barlas et a& based on an evalua- 
tion of graphene's exchange and random phase approx- 
imation (RPA) correlation energies. Moreover, Sheehy 
and Schmalian— by exploiting the proximity to relativis- 
tic electron quantum critical point, derived explicit ex- 
pressions for the temperature and density dependence of 
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the compressibility properties of graphene. All these the- 
oretical efforts have been carried out for clean systems. 
Since disorder is unavoidable in any material, there has 
been great interest in trying to understand how disorder 
affects the physics of electrons in material science spe- 
cially here in graphene and its transport properties. 

Our aim in this work is to study the ground-state prop- 
erties in the presence of electron-impurity and electron- 
electron interactions. For this purpose, we use the self- 
consistent theory of Gotze 22 to calculate the scattering 
rate, ground-state energy and the compressibility of the 
system at the level of RPA including disorder effects. 
Our calculation is in the same spirit of our earlier work 
on conventional 2DEG. 16 We note that recent work of 
Adam et alJ& also use a self-consistent approach where 
the impurity scattering by the charge carriers is treated 
self-consistently in the RPA and the static conductiv- 
ity is calculated in the Boltzmann kinetic theory. Thus, 
the main difference between the present work and that 
of Adam et a?— i is that we are interested in a thermo- 
dynamic quantity (compressibility) whereas the latter is 
aimed at calculating a transport property (conductivity). 
We also remark that direct solution of Dirac equation 
for Dirac-like electrons incorporating the charge impuri- 
ties has been discussed by Novikov22- and the validity of 
the Born approximation is seriously questioned. Similar 
work has been carried out by Pereira et al— in which 
they studied the problem of a Coulomb charge and cal- 
culated the local density of state and local charge by 
solving the Dirac equation. They found new character- 
istics of bound states and strong renormalization of the 
van Hove singularities in the lattice description that are 
beyond the Dirac equation. 

In this work, we consider the charged impurity and the 
surface-roughness potentials which are established exper- 
imentally 2 ^^ to be important. It has been demonstrated 
that a short-range scattering potential is irrelevant for 
electronic properties of graphen o 10 i 26 . We have used the 
same method 16 i 27 to investigate some properties of the 
conventional 2DEG. In this paper, we point out the differ- 
ences between the graphene and conventional 2DEG due 
to disorder effects. The scattering rate behavior within 
our self-consistent theory shows that impurity scattering 
cannot localize the carriers in graphene. The effect of 
disorder on spin susceptibility is similar to that on com- 
pressibility and accordingly we will not show any result 
for spin susceptibility. 

The rest of this paper is organized as follows. In Sec. II, 
we introduce the models for self-consistent calculation of 
impurities effect. We then outline the calculation of com- 
pressibility. Section III contains our numerical calcula- 
tions of ground state properties and comparison of mod- 
els with recent experimental measurements. We conclude 
in Sec. IV with a brief summary. 



II. THEORETICAL MODEL 

We consider a system of 2D Dirac-like electrons inter- 
acting via the Coulomb potential e 2 /er and its Fourier 
transform v q = 2ire 2 / (eg) where e is the background di- 
electric constant. The Dirac electron gas Hamiltonian on 
a graphene sheet is given by 



H = vJ2 & >Q [t 3 ® a- • k] ^ k , Q + v q (n q n_ q - N) 

(1) 

where v = 3ta/2 is the Fermi velocity, t is the tight- 
binding hopping integral, a is the spacing of the honey- 
comb lattice, A is the sample area and N is the total 
number operator. Here r 3 is a Pauli matrix that acts 
on K and K' two-degenerate valleys at which ir and ir* 
bands touch and a 1 and a 2 are Pauli matrices that act 
on graphene's pseudospin degrees of freedom. 

A central quantity in the theoretical formulation of 
the many-body effects in Dirac fermions is the dynam- 
ical polarizability tensor (Q) /i / 0) where [i is 
chemical potential. This is defined through the one- 
body noninteracting Green's functions.— The density- 
density response function \^ (Qj ^> A 4 ) °f t he doped two- 
dimensional Dirac electron model was first consider by 
Strung— as a step toward a theory of collective excitations 
in graphite. The Dirac electron %(°^(q, fi,/x) expression 
has been considered recently by us^ and others— i Imple- 
menting the Green's function (k, w, fi) in the calcu- 
lation a closed form expression for (ij A 1 7^ 0) is 
found.— To describe the properties of Dirac electrons we 
define a dimensionless coupling constant a gr — ge 2 /veh 
where g = g v g s — 4 is the valley and spin degeneracy. 

The effect of disorder is to dampen the charge-density 
fluctuations and results to modify the dynamical polariz- 
ability tensor. Within the relaxation time approximation 
the modified x^^q, ifl, //, T) is given by^i 
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in which the strength of damping is represented by T. To 
include the many-body effects, we consider the density- 
density correlation function within the RPA, 



X (°)(q,»r!, M ,r) 

i-«,x<°>(q J tf*.M.r) 



(3) 



As the short-range disorder is shown^ to have negligi- 
ble effect in the transport properties of graphene, we con- 
sider long-ranged charged impurity scattering and surface 
roughness as the main sources of disorder. The latter 
mechanism also known as ripples comes either from ther- 
mal fluctuations or interaction with the substrate^ 2 - The 
disorder averaged surface roughness (ripples) potential 
(SRP) is modeled as 



(\Usurf{q)n = irA 2 h 2 (2Tre 2 n/eYe 



-q 2 A 2 /4 
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where h and A are parameters describing fluctuations 
in the height and width, respectively. We can use the 
experimental results of Meyer aZ.— who estimate A ~ 
10 nm and h ~ 0.5 nm. It is important to point out 
that there are other models to take into account the 
surface-roughness potential. The effect of bending of 
the graphene sheet has been studied by Kim and Cas- 
tro Neto^. This model has two main effects, firstly the 
decrease of the distance between carbon atoms and sec- 
ondly a rotation of the p z orbitals. Due to bending the 
electrons are subject to a potential which depends on 
the structure of the graphene sheet. Another possible 
model is described by Katsnelson and Geim — consider- 
ing the change of in-plane displacements and out-of-plane 
displacements due to the local curvature of a graphene 
sheet. Consequently, the change of the atomic displace- 
ments results to change in nearest-neighbour hopping pa- 
rameters which is equivalent to the appearance of a ran- 
dom gauge field described by a vector potential. These 
different models need to be implemented in our scheme 
and to be checked numerically to assess their validity in 
comparison to the available measurements. 

The charged disorder potential (CDP) is taken to be 
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in which rij is the density of impurities and d is the set- 
back distance from the graphene sheet. 

We use the mode-coupling approximation introduced 
by Gotze 2 ^ to express the total scattering rate in terms 
of the screened disorder potentials 
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where e(q) = I — WgX^(l) * s the static screen- 
ing function and the relaxation function for elec- 
trons scattering from disorder is given as (po(q,iT) — 

[x (0) (q,tT )M )- x (°)(q)]/ir. 

Since the scattering rate T depends on the relaxation 
function <y9o(q, ir), which itself is determined by the dis- 
order included response function, the above equation 
needs to be solved self-consistently to yield eventually 
the scattering rate as a function of the coupling con- 
stant. Note that at the present level of approximation 
(i.e. RPA) the static dielectric function e(q) does not de- 
pend on r. In the conventional 2DEG correlation effects 
beyond the RPA (through the local-field factor) render 
e(q) also r dependent^ 

The ground state energy is calculated using the cou- 
pling constant integration technique, which has the con- 
tributions E tot = E kin +E X + E C . 

The first-order "exchange" contribution per particle is 
given by 
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FIG. 1: (Color online) The scattering rate T as a function of 
the coupling constant a gr for both the charge-disorder poten- 
tial (CDP) and surface roughness potential (SRP) contribu- 
tions. 
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To evaluate the correlation energy in the RPA, we follow 
a standard strategy for uniform continuum models — 



.RPA 



-t>,x<°>(q,in,^r)]} . 



(8) 



Since x^ (q, f2, ^, T) is linearly proportional to q at 
large q and decrease only like w _1 at large u>, accordingly 
the exchange and correlation energy built by Eqs. (7) 
and (8) is divergent^. In order to improve conver- 
gence, it is convenient at this point to add and subtract 
X^ -* (q, i-i — 0, 2r) inside the frequency integral and 
regularize — the exchange and correlation energy. There- 
fore, these ultraviolet divergences can be cured calculat- 
ing 

I f d 2 a r + °° 
S£ ^^J(2^^l ^X (0) (q,^,r) (9) 
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where is the difference between the doped (/j ^= 0) 

and undoped (fi — 0) polarizability functions. With this 
regularization the q integrals have logarithmic ultraviolet 
divergences^, we can introduce an ultraviolet cutoff for 
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the wave vector integrals k c = Kkp which is the order of 
the inverse lattice spacing and A is dimensionlcss quan- 
tity. Fermi momentum is related to density as given by 
kp = {Aim/ g) 1 / 2 . Once the ground state is obtained the 
compressibility n can easily be calculated from 

K -^ni^^±, (11) 

where the total ground-state energy is given by fetot = 
fckin + Se x + fc^ PA . Here the zeroth-order kinetic contri- 
bution to the ground-state energy is fekm = §£f- we con- 
sider the dimensionless ratio k/kq where kq = 2/(nep) is 
the compressibility of the noninteracting system. 
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FIG. 2: (Color online) (a): The correlation energy 5e c as a 
function of the coupling constant a gr for cut-off value A = 
k c /kF = 50. (b): The exchange energy 8e x as a function of 
the coupling constant a gr for cut-off value A = 50. Results 
of fixed F values are compared to those calculated within the 
mode-coupling approximation. 



III. NUMERICAL RESULTS 

In this section we present our calculations for ground 
state properties of graphene in present of impurities that 
we model them as mentioned above. The inverse com- 
pressibility \/(n 2 K} is calculated by using the theoretical 
models described above and compare them with the re- 
cent experimental measurements. In all numerical calcu- 



lations we consider d = 0.5 nm. Electron density is taken 
to be 1 x 10 12 citT 2 for Figs. 1-3. 

Increasing disorder (increasing n.j or decreasing d 
for charge-disorder potential or increasing h for sur- 
face roughness potential) decrease the x^HQi ^> T) as 
the scattering rate T gets bigger. Thus, decreasing 
X^°Hl, r) ( or increasing correlation effects) results 
in a stronger disorder potential. Despite T increases with 
increasing a gr , apparently it turns to a saturation limit 
and does not diverge. This behavior is quite different 
than what is seen in conventional 2DEG — when the 
many-body effects influence the scattering rate through 
the local-field factor. In the conventional 2DEG sys- 
tem, at a critical level of disorder this nonlinear feedback 
causes T to increase rapidly and diverge, which is taken 
as an indication for the localization of carriers. However, 
in graphene, our calculations show that the T does not di- 
verge therefore impurities cannot localize carriers and we 
have a weakly localized system in the presence of impu- 
rities compatible with experimental observations^. We 
can understand the saturated behavior of T qualitatively 
as follows. In the context of conventional 2DEG, Mott 
argument says that the mean-free path I in a metal could 
not be shorter than the wavelength A. Since I is propor- 
tional to the inverse of V, for large values of T obtained 
in 2DEG, the electron mean- free path decreases and be- 
comes less than or equal to A. At this point we should 
have a metal-insulator phase transition. In the context of 
graphene, on the other hand, Mott's argument is similar 
to that light does not notice any roughness (one source of 
scattering) on a scale shorter than its wavelength. Conse- 
quentely there is a lower limit for the electron's mean-free 
path in graphene and it turns out that we would have a 
maximum (saturation) value for T. 

The issue of localization in graphene has recently at- 
tracted some attention and the chiral nature of elec- 
tron behavior has been discussed in the literature ! 37 ! 38 
Suzuura and Ando 3 ^ claimed that the quantum correc- 
tion to the conductivity in graphene can differ from what 
is observed in normal 2DEG due to the nature of elas- 
tic scattering in graphene possibly changing the sign of 
the localization correction and turn weak localization into 
weak antilocalization for the region when intervalley scat- 
tering time is much larger than the phase coherence time. 
Further consideration of the behavior of the quantum cor- 
rection to the conductivity in graphene^ conclude that 
this behavior is entirely suppressed due to time-reversal 
symmetry breaking of electronic states around each de- 
generate valley. 

We have found through our calculations that T in- 
creases with increasing n,-/n as a function of a gr . Fig- 
ure 1 shows r for various scattering mechanisms. As it 
is clear, CDP is the dominant mechanism for T in our 
model. The effect of SRP is mostly negligible, except at 
large values of the coupling constant. This finding is to 
be contrasted with the the statement of Martin et aZ.— 
that both substrate induced structural distortions (SRP) 
and chemical doping (CDP) are conceivable sources of 
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density fluctuations. We stress that our model calcula- 
tions indicate that at realistic coupling constant values 
(c.f. Fig. 1) only the charged impurity scattering domi- 
nates. 
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FIG. 3: (Color online) The compressibility k/ko scaled by 
that of a noninteracting clean system as a function of the 
coupling constant a gr for cut-off value A = 50. 
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FIG. 4: (Color online) The inverse compressibility [n 2 «:] _1 = 
djj,/dn (in units of meV 10~ 10 cm 2 ) as a function of the elec- 
tron density (in units of 10 12 cm -2 ). The filled squares are 
the experimental data by Martin et aZ.— . 

We have calculated the exchange and correlation ener- 
gies as a function of a gr in the presence of disorder. It is 
found that the disorder effects become more appreciable 
at large coupling constants, within the mode coupling ap- 
proximation. The exchange energy is positive^ because 
our regularization procedure implicitly selects the chem- 
ical potential of undoped graphene as the zero of energy; 
doping either occupies quasiparticle states with positive 
energies or empties quasiparticles with negative energies. 
Figure 2(a) shows the correlation energy 5e c as a func- 
tion of a gr . It appears that the disorder effects become 
more appreciable at large coupling constant. Note that 
5e c has the same density dependence as Se x apart from 



the weak dependence on A. In contrast to the exchange 
energy, Figure 2(b), the correlation energy is negative! 
Figure 3 shows the charge compressibility, k/kq scaled by 
its non-interacting contribution as a function of a gr for 
various models of T. The behavior of k shows some novel 
physics, which is qualitatively different from the physics 
known in the conventional 2DEG. Exchange makes a pos- 
itive contribution to the inverse compressibility and thus 
tends to reduce (rather than enhance) the compressibil- 
ity. On the other hand, correlations make a negative con- 
tribution to the inverse compressibility and thus tends to 
enhance k. In the conventional 2DEG both contribu- 
tions tend to enhance the compressibility. In the case of 
graphene instead, apparently exchange and correlation 
compete with each others in determining the compress- 
ibility of the system. It is interesting to note that similar 
physics is true also in the spin susceptibility 4 . 

In Fig. 4 we compare our theoretical predictions for 
the inverse compressibility of doped graphene with the 
experimental results of Martin et a/^. For deflniteness 
we take A = k c /kp to be such that 7r(Afci?) 2 = (2ir) 2 /Ao, 
where Ao = S\f?iaf j /2 is the area of the unit cell in the 
honeycomb lattice, with ao — 1.42 A the carbon-carbon 
distance. With this choice A ~ (gn -1 \/3/9.09) 1/2 x 10 2 , 
where n is the electron density in units of 10 12 cm -2 . 
Martin et al. ls fitted the experimental inverse compress- 
ibility, (n 2 ^) -1 to the kinetic term using a single param- 
eter Fermi velocity which is larger than the bare Fermi 
velocity. Note that the kinetic term in graphene has the 
same density dependence as the leading exchange and 
correlations terms. 

As it is clear in Fig. 4 the inverse compressibility of 
noninteracting system is below the experimental data. 
By increasing the interaction effects, i.e., increasing the 
coupling constant strength, a gr our theoretical results 
move up. Unfortunately, in the experimental sample, 
the value of a gr is not specified and we considered it to 
be w 1, Therefore, including the exchange-correlation 
effects in our RPA theory, gives results very close to ex- 
perimental data. Furthermore, the results of incorporat- 
ing impurity density, rii — 10 10 cm~ 2 in the system and 
solving the self-consistent equations to obtain the scat- 
tering rate value, yield very good agreement with the 
measured values in the large and middle electron density 
regions. We have examined the inverse compressibility 
by the kinetic term contribution only including a fitting 
value for Fermi velocity and our numerical results are well 
described by a fitting velocity about 1.28 Dp. We would 
stress here that this fitting velocity is different from the 
renormalized velocity defined within the Landau-Fermi 
liquid theory in graphene^ 

In a recent calculation of d[i/dn within the Hartre- 
Fock approximation in graphene where /i is the chemical 
potential and n is the electron density Hwang et al£2- 
stated that correlation and disorder effects would only 
introduce small corrections. This is not, in general, true 
since it has been shown by Barlas et al. 4 that the corre- 
lation effects are essential in the ground-state properties. 
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Although these effects are not significant in very weak 
interaction strength regime and high electron density, in- 
cluding many-body exchange-correlation effects together 
with disorder effect are necessary to get agreement with 
quantities measured in experiments of Martin et alr^ It 
would be useful to carry out further experimental work 
at larger interaction strengths to assess the role played 
by correlation effects. 

IV. CONCLUSION 

We have studied the ground state thermodynamic 
properties of a graphcnc sheet within the random phase 
approximation incorporating the impurities in the sys- 
tem. Our approach is based on a self-consistent calcula- 
tion between impurity effect and many-body electron- 
electron interaction. We have used a model surface 
roughness potential together with the charged disorder 
potential in the system. Our calculations of inverse com- 
pressibility compared with recent experimental results of 



Martin et al} 8 demonstrate the importance of including 
correlation effects together with disorder effects correctly 
in the thermodynamic quantities. 

We remark that in a very small density region, the sys- 
tem is highly inhomogeneous where experimental data 
tends to a constant and the effect of the impurities are 
very essential. A model going beyond the RPA is neces- 
sary to account for increasing correlation effects at low 
density. To describe the experimental data in this region 
more sophisticated theoretical methods which incorpo- 
rate inhomogeneities are needed. One approach would 
be the density-functional theory where Dirac electrons in 
the presence of impurities are considered. 
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